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PHASE  CHARACTERISTICS  AND  TIME  RESPONSES  OF  UNKNOWN  LINEAR  SYSTEMS 
DETERMINED  FROM  MEASURED  CW  AMPLITUDE  DATA 

M.  T.  Ma  and  J.  W.  Adams 

Electromagnetic  Fields  Division 

National  Institute  of  Standards  and  Technology 

Boulder,  CO  80303 

Abstract 

An  alternative  but  simpler  technique  for  calculating  the  complete 
time  and  frequency  characteristics  of  an  unknown  linear  system  from 
the  measured  amplitude  response  to  cw  excitations  is  described.  The 
associated  system  transfer  function  so  determined  may  or  may  not  be  at 
minimum  phase.  A  comparison  of  the  time  responses  shows  the  worst 
case.  Results  also  indicate  that  the  susceptibility  of  the  minimum- 
phase  system  to  damage  by  pulsed  excitation  is  the  greatest  during  the 
initial  period  of  excitation. 

Key  words:  all-pass  function;  Hilbert  transform;  Laplace  transform; 
minimum  phase;  non-minimum  phase;  system  transfer  function;  time 
response . 

1.  INTRODUCTION 

The  time  response  of  a  linear  system  to  a  given  excitation  (whether  it 
be  cw  or  pulse)  can  be  uniquely  determined,  if  the  transfer  function  of  this 
system  is  known,  by  use  of  the  inverse  Laplace  transform.  However,  for 
complicated  systems  involved  in  practical  applications,  the  system  transfer 
function  may  be  unknown.  The  time  response  of  such  a  system  to  pulsed 
excitations  can  be  obtained  only  by  sophisticated  time -domain  measurements 
or  by  derivations  using  measurements  of  the  amplitude  and  phase  responses  to 
cw  excitations.  These  measurement  results  can  then  be  used  to  assess  (a) 
whether  or  not  the  system  may  survive  an  assumed  pulsed  excitation,  and  (b) 
if  not,  what  hardening  is  required  for  the  system  to  withstand  the  threat. 
Such  time -domain  or  frequency- domain  phase  measurements  typically  require 
expensive  equipment  and  special  considerations  on  radiation  hazards  and 
regulatory  compliance  (if  performed  outdoors).  On  the  other  hand,  taking  cw 
amplitude  responses  of  the  same  complicated  system  due  to  low  levels  of 
escitation  is  relatively  easy  and  less  costly.  If  the  measured  cw  amplitude 
could  be  used  to  predict  the  system's  phase  and  thus  its  transfer  function, 
the  time  response  to  a  general  excitation  could  then  be  determined. 

In  this  report,  we  present  an  alternative  but  simpler  technique  to 
determine  the  time  and  frequency  characteristics  of  an  unknown,  linear 
system,  based  only  on  a  given  set  of  cw  amplitude  responses.  The  associated 
system  transfer  function  so  determined  may  or  may  not  be  at  minimum  phase. 
The  time  responses  corresponding  to  the  minimum  phase  and  non-minimum  phases 
thus  obtained  can  then  be  used  to  assess  the  susceptibility  of  the  system  to 
damage  by  pulsed  excitations.  The  theoretical  background  of  relationship 
between  the  amplitude  and  phase  of  a  minimum-phase  transfer  function  is 
outlined  in  section  2.  The  conventional  numerical  approach  practiced  in  the 
past  for  calculating  the  phase  information  and  the  time  response  of  the 
linear  minimum-phase  systems  from  the  measured  cw  amplitude,  and  the 
accuracy   involved   in  this  process  are  reviewed  and  discussed  in  section  3. 


Our  new  approach  for  obtaining  a  set  of  transfer  functions  (both  minimum  and 
non-minimum  phases)  and  the  corresponding  impulse  responses  directly  from 
the  given  cw  amplitudes  is  presented  in  section  4.  This  technique  does  not 
require  knowledge  of  the  necessary  phase  information.  Of  course,  this  phase 
information  is  automatically  obtained  after  the  system's  transfer  function 
is  determined.  It  can  also  be  determined  analytically,  if  it  is  required, 
even  before  the  minimum-phase  transfer  function  is  obtained.  Physical 
meanings  of  the  minimum-phase  results  are  discussed.  A  few  examples  with 
various  orders  of  transfer  functions  and  existing  measured  data  to 
demonstrate  this  idea  are  given.  Energy  contents  associated  with  the  system 
are  discussed  in  section  5.  Some  concluding  remarks  and  extensions  for 
future  work  are  presented  in  the  final  section. 

2.  THEORETICAL  BACKGROUND 

A  well  designed,  stable,  linear  system  with  time- invariant  and  lumped- 
constant  elements  exists  only  when  its  transfer  function  H(s)  has  no  poles 
in  the  right  half  of  the  complex  frequency  s-plane.  In  other  words,  H(s)  is 
analytic  in  Re(s)  >  0,  where  Re  stands  for  the  "real  part  of"  [1].  We 
address  only  the  stable  system  in  this  report,  because  otherwise  the  system 
is  not  useful  in  application.  In  general,  H(s)  is  a  rational  function  of  s 
or  a  ratio  of  two  polynomials  with  the  degree  of  the  numerator  polynomial 
lower  than  the  degree  of  the  denominator  polynomial.  When  this  transfer 
function  is  evaluated  at  the  real  radian  frequency  s  =  jw,  H(jw)  is  a 
complex  function  of  to.  It  then  consists  of  a  real  part  R(w)  and  an 
imaginary  part  X(w)  ,  or  a  magnitude  |H(w)|  and  a  phase  6  (to)  .      That  is, 

H(jco)  =  R(w)  +  j  X(w)  =  |H(jw)|  e"J*(w),  (1) 

where  the  convention  of  assigning  a  minus  sign  to  the  phase  is  adopted.  The 
magnitude  function  |H(jw)|  may  also  be  expressed  in  terms  of  the  attenuation 
function  a(a>)  : 

|H(j»)|  =  e-Q(w).  (2) 

When  H(s)  is  analytic,  the  real  and  imaginary  parts  of  H(jw)  are 
related  by  the  Hilbert  transform  pair  [2], 

X(c)  =  -  **   r  Hr^T  dy, 

and  (3) 

R(«)  =  1   r  f(Y)2  dy. 

In  other  words,  if  one  part  is  given  or  specified,  the  other  part  can  be 
uniquely  determined  by  performing  one  of  the  integrals  shown  in  (3).  The 
complex  function  H(jw)  is  then  completely  obtained.  In  practice,  however, 
(3)  is  not  useful  because  we  cannot  just  measure  R(w)  or  X(w) . 


If,  in  addition,  H(s)  has  no  zeros  in  Re(s)  >  0,  the  transfer  function 
is   said  to  be  at  minimum  phase,   herein  denoted  by  H  (s) .   Under  this 

condition,  the  attenuation  and  phase  functions  are  related  by  another  pair 
of  Hilbert  transforms  [2,  3]: 

,   n  ^nlH  (jy)  I 

-oo  J  -oo   J 

and 

a(W)  =  a(0)  -  f  IT   y(y2(Y).  „*)    dy.  (4b) 

Equation  (4b)  shows  that  the  attenuation  function  can  be  determined 
completely  from  a  given  phase  function  only  when  a(0)  is  also  known.  But, 
for  our  application,  only  (4a)  is  useful  because  we  assume  that  the 
magnitude  (or  attenuation)  function  is  given.  Once  & (w)  is  determined,  the 
entire   complex  H  (jw)   is   obtained   from  (1),  because  |H  (jo) |  or  a(w)  is 

already  known.  The  impulse  response  of  the  system  is  then  traditionally 
calculated  by  the  inverse  Fourier  transform, 

h  (t)  =  7T  r   H  (»  ejWtdt.  (5) 

-oo 

The  system's  time  response  to  a  general  excitation  can  subsequently  be 
computed  by  the  convolution  integral  h  (t)*e(t) ,  where  e(t)  is  an  arbitrary 

excitation,  cw  or  pulse,  applied  to  the  system.  The  success  of  determining 
the  time  response  in  this  case  is  based  on  the  assumption  that  the  system's 
transfer  function  is  at  minimum  phase. 

3.  CONVENTIONAL  APPROACH 

Since  the  type  of  integral  in  (4a)  is  considered  extremely  difficult  to 
calculate,  the  approach  has  been  to  apply  another  transformation  known  as 
the  Wiener-Lee  transform  [2]  to  a(w)  or  |H(ja>)|,  and  then  to  use  numerical 
procedures  [4]  to  obtain  the  necessary  6  (w)  ,  H(ja>)  ,  and  h(t)  . 

When  the  Wiener-Lee  transform  [2], 

w  =  -  tan(5/2),  (6) 

is  applied,  we  are  able  to  transform  the  entire  u>  range,  (-°°,  °°)  into  (-n, 
n)  for  5.  Because  of  this  transformation,  the  original  attenuation  and 
phase   functions,   q(w)   and  9  (u>)  ,  will  become  respectively  a' (8)    and  #'(5). 

Furthermore,   since   the   attenuation   function  a(w)  =  -  in|H(jw)|  =  -  -z   in 

[R2(co)   +  X2(w)]   is   an  even   function,   and   the   phase   function  6  (u>)    = 

tan  [X(w)/R(tc>)  ]  is  an  odd  function,  of  u>  in  (-«,  °°)  ,  their  respective 
transforms  a' (5)  and  0'(8)  are  also  even  and  odd  functions  of  8  in  {-n ,  n) . 
As  such,  they  may  be  expanded  into  Fourier  series, 


a' (8)    =   an  +  a,  cos  8   +  ...  +  a  cos  n5  +  .... 

0  1  n  ' 

and  (7) 

6' (8)    =  b,    sin  5  +  ...  +  b  sin  n5  +  ...  , 

1  n  ' 

where  the  expansion  coefficients  are  determined  by 

a„  ■-■■  -^  J*  a' (8)    dS,  a  =  -  Jn  a' (8)    cos  n8   dS  ,  (8a) 


0    2n  J  n   7r 

-7T  -7T 


and 


b  =  -  Jn   6' (8)    sin  nS  dfi  .  (8b) 


n     7T 

-7T 

When  the   system  under  consideration  is  such  that  h(t)  is  causal  [h(t)  =  0 

when  t  <  0] ,  as  is  usually  so  in  practice,  it  is  known  that  a  and  b   are 

J  r  n      n 

simply  related  by  [2] 

b  =  -  a  .  (9) 

n      n  v  ' 

Therefore,   the  determination  of  a   from  (8a)  automatically  yields  b  ,  which 

in  turn  gives  6' (8),  6 (w) ,  and  thus  H(jw) .  The  impulse  response  is  then 
obtained  from  (5). 

The  justification  for  using  the  Wiener-Lee  transform  and  the  derivation 
as  outlined  above  seems  reasonable  and  straightforward.  The  important 
question  we  should  ask  is  then:  if  the  integral  in  (4a)  is  difficult  to 
compute  before  the  Wiener-Lee  transform  is  applied  because  the  given  |H(jw)| 
or  a(w)  are  very  complicated,  is  it  much  simpler  to  compute  the  Fourier 
series   expansion  coefficients   a  as  required  in  (8a)  after  the  Wiener-Lee 

transform  is  used?  If  so,  what  kind  of  penalty  is  paid  in  terms  of  the 
final  accuracy?  After  all,  the  entire  procedure  requires  (i)  conversion  of 
the  given  a(a>)  data  into  a' (8)  ,    (ii)  numerical  computation  of  a   from  a'  (5)  , 

(iii)  construction  of  a  set  of  6'  (5)  by  including  a  finite  number  of  terms 
in   the   Fourier   sine   series   of   (7)   with  b   =   -   a  ,   (iv)  numerical 

determination  of  H(jw)  based  on  the  given  attenuation  data  and  the  newly 
retrieved  phase   data,  and  (v)  numerical  computation  of  h  (t)  by  performing 

the  discrete  inverse  Fourier  transform.  Each  of  these  five  steps  involves 
approximations  and  thus  introduces  inaccuracies. 

4.  AN  ALTERNATIVE  BUT  SIMPLER  APPROACH 

Here,  we  review  modern  passive  network  theory  together  with  a  re- 
examination of  the  Hilbert  transform  (4a)  to  assess  whether  or  not  a  simple 
and  analytical  method  may  be  devised  to  handle  the  problem  at  hand  without 
going  through  unnecessary  numerical  procedures.  Since  H(s)  for  the  system 
under  study  is  a  rational  function  of  s,  the  square  of  the  magnitude  of 
H(jw) ,  |H(jw)|2,  is  a  ratio  of  two  polynomials  of  even  order  in  w.   As  such, 


the  numerator  and  denominator  polynomials  in  |H(jw)|2  may  contain  one  or 
more,  or  any  combination,  of  the  following  elementary  factors: 

(i)  u>2   +   a2 
and  (10) 

(ii)  w4  +  bw2  +  c2, 

where  a,  b,  and  c  are  real  numbers.  In  addition,  since  |H(ju>)|2  >  0  for  all 
w,  we  require,  for  the  factor  (ii)  in  (10), 

b2  <  4c2.  (11) 

The  number  of  the  above  elementary  factors  contained  in  |H(jw)|2  depends  on 
the  complexity  of  the  transfer  function  being  studied.  The  first-order 
transfer  function  has  just  one  factor  in  form  (i)  in  the  denominator  of 
|H(jw)|2,  and  its  numerator  is  a  constant.  For  a  second-order  transfer 
function,  the  denominator  of  |H(jw)|2  contains  a  factor  in  form  (ii),  or  two 
factors  in  form  (i),  and  the  numerator  has  either  a  factor  in  form  (i)  or  a 
constant.  For  a  third- order  transfer  function,  the  denominator  can  be 
expressed  as  a  product  of  form  (i)  and  form  (ii)  [or  as  a  product  of  three 
factors  in  form  (i)],  and  the  numerator  as  a  factor  of  form  (ii),  form  (i), 
or  a  constant.   Thus,  for  an  nth- order  transfer  function,  the  denominator  of 

|H(jw)|2   has   a   special   form  of  polynomial   such   as   w   +  a   , u   n~"  + 

an./n-4  ♦  .  .  .+  aQ,  and  the  numerator  is  an  even  polynomial  o/the  order 
of  (2n  -  2)  or  lower. 

When  the  magnitude  function  |H(jw)|  can  be  described  by  a  mathematical 
expression,  H(s)  can  be  deduced  directly  and  exactly  from  it  in  a 
straightforward  manner  [5],  without  having  to  find  the  phase  function  first. 
For  example,  consider  the  simplest,  first-order  case  described  by 

|H(>|2  =  ^2  \   a2  ,     a  >  0.  (12) 

We  can  use  the  relationship,  |H(jw)|2  =  H(s)H(-s)|     .   ,  to  obtain 

s  -  jw 

H(s)H(-s)  =  T^L_  =  ___1___  ,  (13) 

yielding 

Hm(s)  =  l/(s  +  a),  (14) 

which  has  a  simple  pole  (s  =  -a)  on  the  left-hand  real  axis.  Since  it  has 
no  zeros ,  the  transfer  function  in  this  particular  form  is  automatically  at 
minimum  phase.  The  associated  impulse  function  is  then  determined  directly 
by   the   inverse  Laplace  transform  as  h  (t)  =  exp(-at),  t  >  0  without  having 

to  find  the  phase  function  first.  It  can  readily  be  verified  that  (14) 
indeed  gives  |H(jw)|2  in  (12).  Note  that  the  other  factor  (-s  +  a)  in  (13) 
is  not  allowable  in  H(s)  because  the  resulting  system  will  then  be  unstable. 


The  expression  in  (12)  may  represent  the  output  of  a  low-pass  filter 
with  its  maximum  at  u  =  0  and  with  its  half -power  (or  3  dB)  point  at  u>  =  a. 
This  same  parameter  "a"  is  also  the  decay  rate  of  the  time  response. 

Clearly,  the  phase  function  associated  with  the  system's  transfer 
function  in  (14)  is, 

0(w)  =  tan_1(w/a) .  (15) 

If,  however,  the  phase  information  is  desired  before  the  system  transfer 
function  in  (14)  is  obtained,  it  can  also  be  found  analytically  by  the 
following  procedure.  Since 

o(w)  =  -  in  |H(jw)|  =  -  \   in  |H(jw)|2  =  J  in(u>2  +  a2),  (16) 

the  required  phase  function  for  the  minimum-phase  case,  in  accordance  with 
(4a) ,  is 


a  r    \        ^    f00     in(y2   +  a2)    dy  _  <±>    r°°  in(y2   +  a2)dy 
H"}    ~   n  J  2(y2    -    w2)  "  n  J  y2    -    w2 


0 


=  tan"V/a),  (17) 

which  agrees  with   (15) .    The   last  step  in  (17)  is  achieved  by  using  the 
definite  integral  [6]: 

rin(p2  +  q2x2)  dx     n    .       -1,  ah  v  ,  .  n  /1QN 

g2x2  -h2 =^htan   <Ji>'   P>^g'h>0-  <18> 

Another  method  for  getting  h  (t)  is  through  the  inverse  Fourier 
integral  given  in  (5),  since  we  already  know  the  complex  H  (jw)  =  |H  (jw) | 
exp[-jfl(w)] . 

That  is , 

h    (t)    =  ^  r  H    (J«)    eJWt   du>  =  r1  T    |H    (J«)  |    J  [wt    "    $  ("}  ]    dc 
m  Z7T  J         m  J  27r  J       '    m  J       ■ 


1    r«o   cos  \cot    -    6  (u>)  1  du> 

'  0 


7T  J  (w2   +   a2)1/2 


1    r«>   a   cos(cc)t)    du>        1    r<»  cj   sin(c^t)    du> 


J.    pco   a   cos  (cot)    dcc>        ip00 

"    7T    i  W2     +    a2  7T    ^ 


u>2    +   a2 


1  -at        1      -at  -at  /ion 

2  e  +  X  e  =  e        »  (iy) 


which  is  the  same  as  obtained  above  by  a  much  simpler  procedure  (inverse 
Laplace  transform).  The  last  line  in  (19)  results  from  the  application  of 
the  following  two  integrals  [6]: 

r   cos (ax)    dx  _  _n      -aB  , 

x2   +  B2        ~   2/3  e  and 


r<»  x   sin(ax)    dx 
J           x2   +  B2 

n      -aB 

"2e 

Re(B)    >  0. 

a  >  0  and  Re(B)    >  0.  (20) 

This  method  of  inverse  Fourier  integral  is  the  one  used  in  the  numerical 
approach  discussed  in  the  last  section. 

We  wish  to  note  another  important  point  about  the  transfer  function. 
Even  though  the  particular  form  in  (14)  is  taken  for  analysis,  it  does  not 
necessarily  imply  that  the  system  under  study  is  minimum  phase.  In  fact, 
the  original  unknown  system  may  possibly  be  represented  by  many  transfer 
functions   with  non-minimum  phases,  because  a  product  of  H  (s)  and  an  all- 

pass   function  yields   a  transfer  function  with  non-minimum  phase  but  still 

with   the   same  I  H(  ico)  1 2 .    That  is,  H  (s)  =  H  (s)  H  (s)  ,  where  H  (s)  is  the 
i  \j  /  i  '   n       ma  n 

non-minimum-phase   transfer   function,   and  H  (s)   is  the  all-pass  function 

a 

defined  as 

H  (s)  =  n  (s  -  a    )/  n(s  +  a  ),  (21) 

with  each  zero  in  the  right-half  s -plane  the  mirror  image  of  the 
corresponding  pole  in  the  left-half  s-plane  [5].  In  (21),  the  symbol  II 
denotes   successive  product,  and  the  parameter  a.    is  either  positive  real  or 

complex.    When  a.  is  complex,  (21)  must  also  have  another  complex  conjugate 

factor.   The  fact  that  |H  (ju>)|2  =  1  explains  why  |H  (ju))|2  =  |H  (jw)|2. 

3.  in  m 

The  corresponding  impulse  response  is  then 

h  (t)  =  h  (t)  *  h  (t),  (22) 

n       m       a 

where  h  (t)   is   the   inverse   Laplace  transform  of  H  (s) ,  and  *  represents 

cl  3. 

convolution  integral.  For  example,  consider  the  case  with  a  =  1  and  H  (s)  = 
(s  -  2)/(s  +2).   We  then  have  h  (t)  =  e"  t,  t  >  0,  and 

h  (t)  =  h  (t)*h  (t)  =  h  (t)*[S(t)  -  4e"2t] 
n       ma       m    L  J 

=  h  (t)  -  r(t),     t  >  0,  (23) 


where 


r(t)  =  4e"t*e"2t  =  4(e-t  -  e"2t).  (24) 


Since   r(t)  >  0  for  all  t  >  0 ,  we  conclude  that  h  (t)  <  h  (t) .   This  implies 

that,  for  this  particular  example,  if  the  minimum-phase  system  can  survive  a 
threat  from  a  given  pulsed  excitation,  the  non-minimum-phase  system  can  also 
successfully  withstand  the  same  threat. 


Before  we  proceed  with  an  example  of  a  second-order  transfer  function, 
we   note   that  the  impulse  response  h  (t)  given  in  (23)  can  also  be  obtained 

by  first  performing  a  partial  fraction  expansion  of  H  (s) , 


n 


H    (s)    =  H    (s)H    (s)    =  t — -   -,  "        -   rr    =      ~~*   .    +  — 7S7  ,                            (25) 

n                 mv    '    av    '         (s   +   l)(s   +2)         s+1        s+2  v/ 

and  then  doing  an  inverse  Laplace  transform  [7]  to  get 

h  (t)  =  -3e_t  +  4e"2t  =  h  (t)  -  r(t) .  (26) 

nv  '                                                m  v  ' 

We  will   now  present   another   example   using  a  second-order  transfer 
function  whose  squared  magnitude  is  given  by 


co2   +   1 
w4  +  w2  +  16  ' 


|H(»I2  =  ,,4  T  ,,2    I  ,c    ,  (27) 


a  plot  of  which  is  shown  in  figure  1  with  its  resonant  frequency  at  to2  =  3 
and  with  a  half -power  bandwidth  of  Au>2  =  12.689.   From  (27),  we  have 

H(s)H(-s)  =  s4  _  s2  -  16  =  (g2  +  3g  +  4)(s2  _  3g  +  4)  ,  (28) 

thus  yielding  the  minimum-phase  transfer  function 

H  (s)  =  2   I  t   TZ  /  >  (29) 

mv  '        sz   +   3s  +  4  '  v  ' 

and  the  simplest  possible  non-minimum-phase  transfer  function 

H  (s)  =  2~1X   1J.  /     ■  (3°) 

nv  /    s-^  +  3s  +  4  v  ' 

As  a  check,  both  (29)  and  (30)  give  the  same  |H(jw)|2  in  (27).  Of  course, 
the  other  possible  non-minimum-phase  transfer  functions  having  the  same 
|H(jw)|2  in  (27)  are  the  product  of  H  (s)  in  (29)  [or  H  (s)  in  (30)]  and  the 

all-pass  function  H  (s)  in  (21). 

cl 

The  phases   associated  with   the  transfer  functions  given  in  (29)  and 
(30)  are  respectively 

0    (w)  =  tan"1[3w/(4  -  w2 ) ]  -  tan"1(w),  (31) 


and 

9    (w)  =  tan"1[3w/(4  -  w2)]  +  tan'^u).  (32) 

The  impulse  responses  of  (29)  and  (30)  can  be  determined,  respectively, 
to  be  [7] 

h  (t)  =  e'L5t(cos  1.32288t  -  0.37796  sin  1.32288t),  (33) 

m 


and 


h  (t)  =  e~1,5t(-  cos  1.32288t  +  1.88982  sin  1.32288t) 
n 

=  hm(t)  -  r(t),  (34) 

where 

r(t)  =  2e"1-5t(cos  1.32288t  -  1.13389  sin  1.32288t)  >  0 

in   0  <   t  <  0.54634.    This   means  that  h  (t)  <  h  (t)  shortly  after  an 

nv  '  mv  '  J 

excitation  is  applied  to  the  system.  The  minimum-phase  response  will  not 
underestimate  the  system  behavior. 

Thus,  we  have  shown  once  again  that  (i)  the  minimum -phase  transfer  function 
for  an  unknown  system  and  the  resultant  impulse  response  can  be  determined 
directly  from  the  given  |H(jw)|2  without  having  to  obtain  the  associated 
phase  function,  and  (ii)  the  minimum-phase  time  response  h  (t)  so  determined 

is  the  most  conservative  estimate  as  far  as  the  initial  threat  to  the  system 
by  an  excitation  is  concerned. 

The  phase  function  corresponding  to  the  minimum  phase  can  also  be 
determined  from  (4a)  before  its  transfer  function  in  (29)  is  obtained.  We 
begin  with  the  denominator  of  (27) , 

w4  +  w2  +  16  =  (w2  +  af)(w2  +  a2,)  ,  (35) 

where  a\  and  a\  are  a  complex  conjugate  pair.  Even  though  it  is 
straightforward  to  determine  a\  and  a2,  [in  this  case,  a\  =  (1  +  jj63)/2,  and 
a2  =  (1  -  j J  63)/2 ] ,  we  only  need  to  know,  as  shown  later,  axa2  and  ax  +  a2 . 
We  know  from  (35)  that  a\  +  a\  =  1  and  a\a\  =  16.  We  obtain  axa2  =  4,  and 
(ax  +  a2)2  =  a\   +  a%   +   2axa2  =1+8=9,  ora1+a2=3.   Then  we  have 

in|H(jw)|    =    -   |   in(w2    +   a?)    -   |  in(w2    +   a|)    +  |   in(w2   +   1).  (36) 

Applying  (4a)  we  obtain 

0(w)  =  61(o))    +    62(u)    -   83(u>),  (37) 


where,    with   the   aid  of    (18), 

in(y2   +  a\)  . 

M«>)    =  *Z  r  —^2 ^2—  dy  =  tan"    (w/ai),  (38) 

n   Q  y  w 

in(y2   +   a2,)  , 

M")    =  ^  r        y2    .    ^2        dy  =   tan"    (»/a2),  (39) 


and 


a    /    n         w    f°°   ln(v2   +   1)  ^_      -1.    .  ..... 

M")    =  ~  J  y2    .    w2      dy  =   tan      (w) .  (40) 


0 
However , 


!  (w)    +   #2(w)    =   tan      (u)/a1)    +   tan      (w/a2) 


tan      [ (a2    +   a2)w/(a1a2    -    o>2 ) ] 


=   tan"1[3w/(4    -    u>2)]  .  (41) 


Thus ,  we  have 


0(w)  =  tan"1[3w/(4  -  w2)]  -  tan'V),  (42) 

which  checks  with  (31) . 

The  impulse  response  given  in  (33)  can  also  be  obtained  by  the  inverse 
Fourier  integral  method,  as  we  demonstrated  in  the  first  example  for  the 
first-order  transfer  function,  with  a  much  more  involved  procedure. 

We  make  two  comments  about  the  examples  presented  above.  First,  in 
obtaining  (38)  and  (39),  we  have  extended  the  application  of  (18)  to  the 
case  with  a  complex  number  for  p.  Evidently,  (18)  is  still  valid  as  long  as 
a  complex  conjugate  pair  appears  together.  This  is  certainly  satisfied  in 
our  study  since  the  given  |H(jw)|2  is  a  rational  function  of  w2 . 

Second,  the  variable  o>  has  been  used  loosely  as  frequency.  In  (27)  we 
have  identified  u>2  =  3  (or  u>  =  J  3)  as  the  resonant  frequency  for  the  second- 
order  case.  Strictly  speaking,  u>  should  be  treated  as  the  normalized  radian 
frequency.  It  can  be  translated  into  any  frequency  of  interest  by  a  simple 
frequency  transformation  [5], 

w'  =  Aw,  (43) 

where  A  is  a  normalization  constant. 
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If,  for  example,  the  actual  resonant  frequency  occurs  at  J  3  MHz,  we 
then  use  A  =  2tt(10)6  to  translate  u>  =  J  3  to  u>'  =  J  3  x  27r(10)6.  Substituting 
u>  =  u>' /A  into  (27)  gives  a  new  squared  magnitude  function, 


A2(cj'2  +  A2)      ,  (44) 


IgO')|2  -  ,y4  +  A2u>'2  +  16A4 


yielding  the  |G(jw')|2  vs  w'2  curve  identical  to  that  shown  in  figure  1. 
The   corresponding  minimum-phase   transfer   function,  impulse  response,  and 
phase  function  then  become  respectively: 


A(s'  +  A) 
nT   '    s' 


G-(s')  -.■?';  3L-+  4A»  •  <45) 


-1  SAr 
g  (t)  =  A  e   •    (cos  1.32288At  -  0.37796  sin  1.32288At),        (46) 


and 

0(w')  =  tan'^Aw'/OA2  .  w'2)]  -  tan'V'  /A)  .  (47) 

The  same  principle  and  procedures  apply  to  higher-order  cases.  The 
results   on  extracting  H  (s) ,  h  (t) ,  and  6 (w)  from  a  given  |H(jw)|2  and  the 

procedures  involved,  as  demonstrated  in  the  two  examples,  are  exact,  and  no 
approximations  have  been  exercised.  Approximations  are  needed  only  when  an 
analytical  expression  for  |H(jw)|2  is  unknown.  In  practice,  when  only  the 
measured  data  of  |H(jw)|  are  available,  all  we  have  to  do  is  to  obtain  an 
experimental  curve  for  |H(jw)|2  and  then  to  approximate  it  with  a 
realizable,  mathematical  expression.  By  realizability ,  we  mean  (i)  |H(jw)|2 
=  N(c<;2)/D(w2  )  >  0  for  all  w,  (ii)  the  numerator  and  denominator  polynomials, 
N(w2)  and  D(w2)  are  in  even  orders  of  to  with  some  restrictions  on  their  real 
coefficients,  and  (iii)  the  degree  of  N(w2)  is  lower  than  the  degree  of 
D(w2)  by  at  least  2.  This  can  be  achieved  by  the  available  approximation 
theory  and  curve -fitting  techniques.  In  fact,  this  is  the  only  step 
involving  approximations.   The  remaining  procedures  will  be  exact. 

The  last  example  to  demonstrate  the  applicability  of  the  approach 
presented  in  this  section  is  now  given,  and  is  based  on  the  measured  data 
shown  in  figure  2  (solid  curve) .  The  data  represent  normalized,  vertically 
polarized  electric  fields  reflected  from  a  helicopter  when  it  is  irradiated 
by  an  impulse.  By  examining  the  solid  curve  in  figure  2  more  carefully,  we 
notice  at  least  three  resonant  frequencies  at  16.5  MHz,  25.6  MHz,  and  54.3 
MHz.  Obviously,  this  system  is  too  complicated  to  be  represented  by  a 
simple  second-order  transfer  function.  However,  for  demonstration  purpose, 
we  only  consider  the  most  important  resonant  frequency  at  25.6  MHz.  This 
means  that  we  treat  the  dashed  curve  in  figure  2  as  "the  given"  one,  and 
then  try  to  approximate  it  by  a  second-order  transfer  function  with  minimum 
phase  in  the  form  of 

H  (s)  =  B/(s2  +  as  +  b),  (48) 

where  B,  a,  and  b  are  three  unknown  constants  to  be  determined. 
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From  (48) ,  we  have 

"Hm(Jw)l2  =  w*  -  (2b  ^a2)  u>2  +  b2  •  <49> 

which  has  the  resonant  frequency  at  u>%   =   b  -  a2/2. 

Expressing  (49)  in  terms  of  u>0    gives 

lHm(»l2  ~u<    -  ill   a,2  +b2  '  <50> 

Two  essential  conditions  are  then  imposed  on  (50) : 

<4  =  b  -  a2/2  =  (2tt  x  25.6  x  106)2  ,  (51) 

and 

|H  (jw0)l2  =  532    (°r  34.48  dB  as  given  in  fig.  2)  (52) 

A  third  condition  may  be  imposed  at  either  (i)  20  MHz  or  (ii)  30  MHz. 
Indeed,  any  other  condition  based  on  a  physical  reasoning  may  be  used.  If 
we  choose  (i)  with  wa  =  2n   x  20  x  (10)6  =  1.2566  x  108 ,  and  use 

|Hm(jWi)|2  =  132    (or  22.28  dB  as  indicated  in  fig.  2),        (53) 

we  have  the  solutions  from  (51) ,  (52) ,  and  (53) : 

B  =  13.5185  x  1016  ,    a  =  0.1584  x  108 ,    b  =  2.5998  x  1016.    (54) 

The   |H  (jw)|2   in   (50)   with  the  values   of  B,  a,  and  b  given  in  (54)  is 

plotted  in  figure  2  as  the  curve  marked  with  dots.  Not  surprisingly,  the 
frequency  region  from  10  MHz  to  25.6  MHz  is  very  accurately  matched  to  the 
dashed  curve,  while  the  frequency  region  greater  than  25.6  MHz  is  not. 

If  we  choose  (ii)  for  matching  at  30  MHz  (which  is  a  3-dB  point),  we  have 
w2  =  1.8850  x  108.   The  condition 

|Hm(j^2) |2  =  532/2    (or  31.48  dB),  (55) 

with  (51)  and  (52),  yields 

B  =  51.1874  x  1016,    a  =  0.5906  x  108 ,    b  =  2.7616  x  1016.     (56) 

The   |H  (jw)|2   in  (50)  with  this  new  set  of  B,  a,  and  b  values  is  presented 

in  figure  2  as  the  curve  marked  with  crosses.  Here  we  see  a  better 
approximation  for  frequencies  higher  than  the  resonant  frequency,  but  a 
poorer  approximation  for  lower  frequencies. 
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Even  though  these  two  approximations  may  not  be  wholly  satisfactory, 
they  are  still  considered  effective  because  the  complicated  experimental 
curve  is  approximated  by  a  simple  second-order  transfer  function.  A  higher- 
order  (4th  or  6th)  transfer  function,  or  even  a  different  kind  of  a  second- 
order   transfer   function  such  as   H  (s)   =  B(s  +  c)/(s2  +  as  +  b)  with  an 

additional  control  parameter  c  will  certainly  give  a  better  approximation. 
The  use  of  higher-order  transfer  functions  and  additional  control  parameters 
should  be  topics  for  future  work.  The  objective  of  this  report  is  primarily 
to  demonstrate  the  feasibility  of  solving  this  kind  of  problem  with  a 
simpler  approach. 

The  impulse  responses  for  the  two  approximations  under  discussion  are: 

(i)  with  the  matching  point  at  20  MHz, 

h  ,(t)  =  8.4347  x  108  e"  ait  sin  ^t,     t  >  0,  (57) 


ml 


with  ai  =  0.0792  x  108 ,  and  fi1    =  1.6104  x  108 


and 

(ii)  with  the  matching  point  at  30  MHz, 

hm2(t)  =  31.0755  x  108  e"  a2<:  sin  /32t  ,  t  >  0  ,  (58) 

with  a2   =   0.2953  x  108 ,  and  02   =  1.6354  x  108  . 

The  normalized  time  responses  are  shown  in  figure  3  with  curves  marked 
respectively  with  dots  and  crosses,  along  with  the  normalized  measured  time 
response  (the  unmarked  curve).  On  the  first  look,  the  predicted  minimum- 
phase  time  responses  clearly  do  not  agree  well  with  the  measured  curve. 
This  implies  that  the  actual  system  (helicopter  plus  its  attached 
instruments  and  equipment)  may  have  a  non- minimum -phase  transfer  function 
with  one  or  more  zeros  in  the  right  half  s -plane.  This  fact  is  also 
indicated  by  the  initial  negative  response  appearing  in  the  experimental 
curve.    On   the   other  hand,  the  general  shape  of  a  damped  sinusoid  with  a 

_  Q  _  O 

period  between  3 . 84  x  (10)  and  3.90  x  (10)  seconds  agrees  well  with  the 
shape  of  the  experimental  curve.   In  addition,  the  rate  of  decay  for  h  «(t) 

with  a2      =  0.2953   x   108   also   approaches   the  experimental  curve.   This 

indicates   that  matching  at   30  MHz  not  only  improves  |H  (jw)|2  at  higher 

frequencies,  but  also  yields  a  better  approximation  to  the  measured  time 
response.  This  is  desirable  because  the  receiving  antenna  used  in  the 
actual  measurement  also  works  more  accurately  in  the  upper  frequency  range . 

The  corresponding  phases  are: 

e1(u>)    =   tan"1[0.1584  x   108   w/(2.5998   x   1016    -    u>2 )  ]  ,  (59) 

and 
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B2(u>)    =   tan_1[0.5906   x   108   w/(2.7616   x   1016    -    u>2  )  ]  .  (60) 

5.  CONSIDERATION  OF  ENERGY  CONTENTS 

To  assess  the  ability  of  a  system  to  withstand  damage  from  an 
excitation,  it  is  often  useful  to  compute  the  energy  associated  with  an 
impulse  response.  Indeed,  if  h(t)  represents  a  voltage  waveform  across  a  1 
Q   resistor,  the  quantity 

E  =  J°°h2(t)  dt  (61) 

0 

equals   the   total   energy  delivered  to  the  resistor  by  the  excitation  [2] . 
Equation  (61)  also  represents  the  area  under  the  curve  h2(t). 

The  energy  E  may  also  be  computed,  in  view  of  Parseval's  theorem  [2],  by 
E  =  ^  r    |H(j")|2  dw.  (62) 

-00 

Since   the  minimum-phase   impulse   response  h  (t)   and  the  associated  non- 
minimum-phase   impulse   response  h  (t)  =  h  (t)*h  (t)  have  the  same  |H(jw)|2, 

XT  m         a. 

their  respective  total  energies  [in  0  <  t  <  °°]  are  equal  even  though  h  (t)  < 
h  (t)  during  the  initial  period  near  t  =  0   as  demonstrated  in  section  4. 

For  the  first  example  considered  in  section  4,  we  have 

IH    (jw)|2    =   l/(w2   +  a2),  and  h    (t)    =   e"at. 

i    m\j    /i  /  /  m 

The  energy  is  then,  from  (62), 

E  -  ^  r  -a^-2  =  i  r  -TW  =  V(2a)  ;  (63) 

-00  0 


or,  from  (61) , 


E  =  /VCt)    dt  =  f°e-Zat   dt  =   l/(2a).  (64) 

0      m  0 


Referring  to  the  example  with  a  non-minimum-phase  impulse  response  presented 

in  (23),  we  have  h  (t)  =  e_t  (with  a  =  1) ,  h  (t)  =  6(t)  -  4e"2t,  and  h  (t)  = 
\   /  .  m  a  n 

-t      -2t 
h  (t)  *  h  (t)  =  -3e   +  4e    .   The  energy  is  then 
mv  '  av  bJ 

2t    nl     -3t    ,,  -4t, 


E  =  JN°°h2(t)  dt  =  f°(9e'Zt    -    24e"-3t  +  16e   C)  dt 


0   n         0 
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=  1/2  =  l/(2a)  with  a  =  1,  (65) 

which  is  indeed  the  same  as  in  (63)  and  (64) . 

For  the  second  example  as  given  in  (27) ,  we  have 

do  a*dw 

„   _  _1  f<o  (u>2    +    l)do)      1  r«>      1 1  r<»     1 

=  2,  Jm  w*  +  u>*   +  16  "  *  J   2(«»  +  ^)     *  J  2(w2  +  ^ 

=  |   {aj/J^!   +  at/J/8t)    =   5/24,  (66) 

where 

a,    -  LLi/lM    _  and  /Ji   _  Lij!63  (6?) 

The   same   result   can  be   obtained  by   computing 

E  =  fV(t)    dt  or  E  =  fV(t)    dt,  (68) 

0  0 

where  h  (t)  and  h  (t)  are  given  respectively  in  (33)  and  (34). 

For  our  third  example  involving  experimental  data,  the  energies  are: 
(i)  when  the  frequency  point  is  matched  at  20  MHz, 

E,  =  f°0h21(t)  dt  =  2.229  x  (10)10, 
1   £  ml   '  \   /   . 

and  (ii)  when  the  frequency  point  is  matched  at  30  MHz, 

E2  =  rV2(t)  dt  =  7.978  x  (10)10, 

where  h  -,  (t)  and  h  ?(t)  are  given  respectively  in  (57)  and  (58). 

6.  CONCLUSIONS 

We  have  shown  an  alternative  and  simpler  method  to  determine  the 
complete  transfer  function  (including  the  associated  phase  information)  and 
the  impulse  response  for  an  unknown  linear  system  from  a  given  cw  amplitude 
response,  without  having  to  retrieve  the  phase  function  first.  The 
solutions  (both  minimum  and  non-minimum  phases)  are  exact  only  if  an 
analytical  expression  for  the  cw  amplitude  response  is  known.  If  not, 
approximations  can  be  used  to  derive  a  realizable  squared  magnitude 
describing  the  cw  system  response.  The  remaining  procedures  for  obtaining 
the  complete  system  characteristics  are  exact.  Three  examples,  using  known 
analytical  expressions  and  measured  data  for  the  cw  magnitude,  have  been 
examined  and  shown  to  validate  the  proposed  technique.  We  also  have  shown 
that   the  minimum-phase   time   response   represents   the  most  conservative 
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estimate  as  far  as  the  response  to  an  initial  threat  to  the  system  by  an 
excitation  is  concerned.  Examination  of  the  basic  properties  of  transfer 
functions  with  higher  orders  and  with  different  forms  for  each  order,  and 
study  of  accuracies  involved  in  the  approximation  (when  the  measured  data 
are  available)  to  cover  more  complicated  cases  should  be  topics  for  future 
research. 
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